home *** CD-ROM | disk | FTP | other *** search
/ Complete Linux / Complete Linux.iso / docs / apps / circuits / spice2g6.z / spice2g6 / spice / Fortran / scale.f < prev    next >
Encoding:
Text File  |  1989-02-03  |  2.0 KB  |  68 lines

  1.       subroutine scale(xmin,xmax,n,xminp,xmaxp,del)
  2.       implicit double precision (a-h,o-z)
  3. c
  4. c     this routine determines the 'optimal' scale to use for the plot of
  5. c some output variable.
  6. c
  7. c
  8. c  adapted from algorithm 463 of 'collected algorithms of the cacm'
  9. c
  10. c spice version 2g.6  sccsid=knstnt 3/15/83
  11.       common /knstnt/ twopi,xlog2,xlog10,root2,rad,boltz,charge,ctok,
  12.      1   gmin,reltol,abstol,vntol,trtol,chgtol,eps0,epssil,epsox,
  13.      2   pivtol,pivrel
  14.       integer xxor
  15.       dimension vint(5)
  16.       data vint / 1.0d0,2.0d0,5.0d0,10.0d0,20.0d0 /
  17.       data eps / 1.0d-12 /
  18. c
  19. c
  20. c...  trap too-small data spread
  21. c***********************************************************
  22. c  temporily check 'equality' this way
  23.       if(xmin.eq.0.0d0.and.xmax.eq.0.0d0) go to 4
  24.       if(dabs((xmax-xmin)/dmax1(dabs(xmin),dabs(xmax))).ge.1.0d-4)
  25.      1  go to 10
  26.     4 continue
  27.       if (xmin.ge.0.0d0) go to 5
  28.       xmax=0.5d0*xmin+eps
  29.       xmin=1.5d0*xmin-eps
  30.       go to 10
  31.     5 xmax=1.5d0*xmin+eps
  32.       xmin=0.5d0*xmin-eps
  33. c...  find approximate interval size, normalized to [1,10]
  34.    10 a=(xmax-xmin)/dble(n)
  35.       nal=idint(dlog10(a))
  36.       if (a.lt.1.0d0) nal=nal-1
  37.       xfact=dexp(xlog10*dble(nal))
  38.       b=a/xfact
  39. c...  find closest permissible interval size
  40.       do 20 i=1,3
  41.       if (b.lt.(vint(i)+eps)) go to 30
  42.    20 continue
  43.       i=4
  44. c...  compute interval size
  45.    30 del=vint(i)*xfact
  46.       fm1=xmin/del
  47.       m1=fm1
  48.       if (fm1.lt.0.0d0) m1=m1-1
  49.       if (dabs(dble(m1)+1.0d0-fm1).lt.eps) m1=m1+1
  50. c...  compute new maximum and minimum limits
  51.       xminp=del*dble(m1)
  52.       fm2=xmax/del
  53.       m2=fm2+1.0d0
  54.       if (fm2.lt.(-1.0d0)) m2=m2-1
  55.       if (dabs(fm2+1.0d0-dble(m2)).lt.eps) m2=m2-1
  56.       xmaxp=del*dble(m2)
  57.       np=m2-m1
  58. c...  check whether another loop required
  59.       if (np.le.n) go to 40
  60.       i=i+1
  61.       go to 30
  62. c...  do final adjustments and correct for roundoff error(s)
  63.    40 nx=(n-np)/2
  64.       xminp=dmin1(xmin,xminp-dble(nx)*del)
  65.       xmaxp=dmax1(xmax,xminp+dble(n)*del)
  66.       return
  67.       end
  68.